Matlab Codes

Chapter 3

Matlab code 3.1: Matlab file “Figure 3-6.m”

%--------------------------------------------------------------------

% This code can be used to generate Figure 3.6

%--------------------------------------------------------------------

 

clear all;

close all;

%% the theoretical curve of the rotating target micro-Doppler characteristic(LFM radar)

c = 3e8;

j = sqrt(-1);

fc = 10e9; % carrier frequency of transmitted signal

B = 200e6; % bandwidth

cord = 1000*[3 4 5]; % coordinates of local coordinate system's origin in the radar coordinate system

colo = [0 0 0]; % coordinates of radar in the radar coordinate system

R0 = cord-colo;

n = R0/sqrt(sum(R0.^2)); % unit vector of Line-of-Sight(LOS)

ae = [0 pi/4 pi/5]; % Euler angle

tgt = [0 0 0;3 1.5 1.5;-3 -1.5 -1.5]; % three point scatterers

w = [pi 2*pi pi]; % angular velocity of rotation

T = 0.8165; % rotating period

ri1 = [cos(ae(3)) sin(ae(3)) 0;-sin(ae(3)) cos(ae(3)) 0;0 0 1];

ri2 = [1 0 0;0 cos(ae(2)) sin(ae(2));0 -sin(ae(2)) cos(ae(2))];

ri3 = [cos(ae(1)) sin(ae(1)) 0;-sin(ae(1)) cos(ae(1)) 0;0 0 1];

ri = ri1*ri2*ri3; % initial rotation matrix

tn = length(w); % number of scatterers

r0 = zeros(3);

for i = 1:tn

    r0(i,:) = tgt(i,:)-colo; % scatterers in the local coordinate system

end

r0r = ri*r0'; % scatterers in the reference coordinate system

w = ri*w'; % angular velocity of rotation in the reference coordinate system

omega = sqrt(sum(w.^2));

we = w/omega; % unit vector of rotation

wr = [0 -we(3) we(2);we(3) 0 -we(1);-we(2) we(1) 0]; % skew symmetric matrix

v = 0; % translational velocity of target

t = 3; % radar illumimated time

tp = 5e-6; % pulse duration

prf = 1000; % pulse repetition frequency

pri = 1/prf; % pulse repetition interval

dt = 0:pri:t-pri; % slow time sampling interval

m = length(dt);

rdt = zeros(tn,m);

for i = 1:m

    rdt(:,i) = ((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r)'*n';

end

plot(dt,rdt)

xlabel('Slow time (s)')

ylabel('Range (m)')

axis([0,3,-6,6])

 

%% the simulation result of the rotating target micro-Doppler characteristic(LFM radar)

xmin = sqrt(sum((cord-colo).^2)); % the minmum detectable distance

xmax = sqrt(sum((cord-colo).^2))+sqrt(sum((tgt(2,:)-colo).^2)); % the maximum detectable distance

ts = xmin/c; % the starting point of the time sampling

te = xmax/c+tp; % the ending point of the time sampling

mu = B/tp; % frequency modulation rate

fs = B; % fast time sampling frequency

kdt = 1/fs; % fast time sampling interval

nt = 2*ceil((te-ts)/(2*kdt)); % number of fast time sampling

df = fs/nt; % number of frequency sampling

Rr = zeros(length(tgt(1,:)),length(dt)); % distance between the scatterers and radar

s = zeros(nt,length(dt)); % echo signal matrix

for i = 1:m

    for n = 1:length(tgt(:,1))

        Rr(n,i) = sqrt((R0'+v*dt(i)+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n)))'...

                      *(R0'+v*dt(i)+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n))));       

    end

end

tk = ts+kdt*(0:nt-1); % fast time

for n = 1:length(tgt(:,1))

    td = tk(:)*ones(1,m)-2*ones(nt,1)*Rr(n,:)/c; % echo delay

    s = s+exp(j*2*pi*fc*td+j*pi*mu*td.^2);

end

Rr0 = sqrt(sum(R0.^2));

td0 = tk(:)*ones(1,m)-2*Rr0/c;

s0 = exp(j*2*pi*fc*td0+j*pi*mu*td0.^2); % reference signal

sdt = s.*conj(s0); % dechirp

sdf = fftshift(fft(sdt),1);

x = dt;

y = -(df*(-nt/2:nt/2-1))*c/(2*mu); % turn the fast time frequency-slow time plane to the range-slow time plane

figure

contour(x,y,abs(sdf),10)

xlabel('Slow time (s)')

ylabel('Range (m)')

axis([0,3,-6,6])